#include <iostream>
#include "utils.h"

extern float XMIN, YMIN, ZMIN;
extern float XMAX, YMAX, ZMAX;
extern int DX, DY, DZ;
extern int DXY, DXYZ;
extern unsigned int NUMBINS;
extern float MAXPROBE;
extern float GRID;
extern float GRIDVOL;
extern float WATER_RES;
extern float CUTOFF;
extern char XYZRFILE[256];

int main(int argc, char *argv[]) {
  cerr << endl;

  COMPILE_INFO;
  CITATION;

// ****************************************************
// INITIALIZATION
// ****************************************************

//HEADER INFO
  char file[256]; file[0] = '\0';
  char ezdfile[256]; ezdfile[0] = '\0';
  char pdbfile[256]; pdbfile[0] = '\0';
  double BIGPROBE=10.0;
  double SMPROBE=1.5;
  double TRIMPROBE=1.5; //HEADER INFO

  while(argc > 1 && argv[1][0] == '-') {
    if(argv[1][1] == 'i') {
      sprintf(file,&argv[2][0]);
    } else if(argv[1][1] == 's') {
      SMPROBE = atof(&argv[2][0]);
    } else if(argv[1][1] == 'b') {
      BIGPROBE = atof(&argv[2][0]);
    } else if(argv[1][1] == 't') {
      TRIMPROBE = atof(&argv[2][0]);
    } else if(argv[1][1] == 'g') {
      GRID = atof(&argv[2][0]);
    } else if(argv[1][1] == 'o') {
      sprintf(pdbfile,&argv[2][0]);
    } else if(argv[1][1] == 'e') {
      sprintf(ezdfile,&argv[2][0]);
    } else if(argv[1][1] == 'h') {
      cerr << "./FsvCalc.exe -i <file> -b <big_probe> -s <small_probe> " << endl
	<< "\t-t <trim probe> -e <EZD outfile> -o <PDB outfile> -g <gridspace>" << endl;
      cerr << "FsvCalc.exe -- Calculates the Fractional Solvent Volume" << endl;
      cerr << endl;
      return 1;
    }
    --argc; --argc;
    ++argv; ++argv;
  }


//INITIALIZE GRID
  finalGridDims(BIGPROBE);

//HEADER CHECK
  cerr << "Grid Spacing: " << GRID << endl;
  cerr << "Resolution:      " << int(1000.0/float(GRIDVOL))/1000.0 << " voxels per A^3" << endl;
  cerr << "Resolution:      " << int(11494.0/float(GRIDVOL))/1000.0 << " voxels per water molecule" << endl;
  cerr << "Complexity:      " << int(8000000/float(GRIDVOL))/1000.0 << endl;
  cerr << "Input file:   " << file << endl;

//FIRST PASS, MINMAX
  int numatoms = read_NumAtoms(file);

//CHECK LIMITS & SIZE
  assignLimits();

// ****************************************************
// STARTING FIRST FILE
// ****************************************************


//GET SHELL
  gridpt *shell;
  shell = (gridpt*) malloc (NUMBINS);
  if (shell==NULL) { cerr << "GRID IS NULL" << endl; exit (1); }
  zeroGrid(shell);
  int shellvol = get_ExcludeGrid_fromFile(numatoms,BIGPROBE,file,shell);

//INIT NEW chanACC GRID
  cerr << "Trimming Radius: " << TRIMPROBE << endl;
  gridpt *smShell;
  smShell = (gridpt*) malloc (NUMBINS);
  if (smShell==NULL) { cerr << "GRID IS NULL" << endl; exit (1); }

//COPY AND TRUNCATE (IF NECESSARY)
  copyGrid(shell,smShell);
  if(TRIMPROBE > 0) {
    trun_ExcludeGrid(TRIMPROBE,shell,smShell);
  }
  free (shell);

//COPY SMSHELL INTO CHANACC
  gridpt *chanACC;
  chanACC = (gridpt*) malloc (NUMBINS);
  if (chanACC==NULL) { cerr << "GRID IS NULL" << endl; exit (1); }
  int smshellvol = copyGrid(smShell,chanACC);

//SUBTRACT WATER_ACC FROM SHELL TO GET ACC CHANNELS
//GET ACCESS WATER VOLUME
  gridpt *waterACC;
  waterACC = (gridpt*) malloc (NUMBINS);
  if (waterACC==NULL) { cerr << "GRID IS NULL" << endl; exit (1); }
  zeroGrid(waterACC);
  fill_AccessGrid_fromFile(numatoms,SMPROBE,file,waterACC);
  subt_Grids(chanACC,waterACC);
  free (waterACC);

//INIT NEW chanEXC GRID
  gridpt *chanEXC;
  chanEXC = (gridpt*) malloc (NUMBINS);
  if (chanEXC==NULL) { cerr << "GRID IS NULL" << endl; exit (1); }

//GROW EXCLUDED SURFACE FROM ACCESSIBLE
  copyGrid(chanACC,chanEXC);
  grow_ExcludeGrid(SMPROBE,chanACC,chanEXC);
  free (chanACC);

//INTERSECT
  intersect_Grids(chanEXC,smShell); //modifies chanACC
  free (smShell);

//OUTPUT
  int chanvol = countGrid(chanEXC);
  cout << "BIG PROBE (R = " << BIGPROBE << ") VOLUME: ";
  printVolCout(shellvol);
  cout << endl << "TRIMMED (BY R = " << TRIMPROBE << ") VOLUME: ";
  printVolCout(smshellvol);
  cout << endl << "SOLVENT (R = " << SMPROBE << ") VOLUME: ";
  printVolCout(chanvol);
  cout << endl << "FRACTIONAL SOLVENT VOLUME: " << double(100.0*chanvol) / double(shellvol) << "%" << endl << endl;

  if(pdbfile[0] != '\0') {
    write_SurfPDB(chanEXC,pdbfile);
  }
  if(ezdfile[0] != '\0') {
    write_HalfEZD(chanEXC,ezdfile);
  }
  free (chanEXC);

  cerr << endl << "Program Completed Sucessfully" << endl << endl;
  return 0;
};

